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Abstract. 


This  paper  presents  algorithms  for  the  assignment  problem,  the  transportation  problem  and  the 
minimum  cost  flow  problem  of  operations  research.  The  algorithms  find  a  minimum  cost  solution, 
but  run  in  time  close  to  the  best-known  bounds  for  the  corresponding  problems  without  costs.  For 
example,  the  assignment  problem  (equivalently,  minimum  cost  matching  on  a  bipartite  graph)  can 
be  solved  in  0{s/nm  log(nJV))  time,  where  n,  m  and  N  denote  the  number  of  vertices,  number 
of  edges  and  largest  magnitude  of  a  cost;  costs  are  assumed  to  be  integral.  The  algorithms  work 
by  scaling.  As  in  the  work  of  Goldberg  and  Tarjan  [Go,  GoT87a-b],  in  each  scaled  problem  an 
approximate  optimum  solution  is  found,  rather  than  an  exact  optimum. 
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1.  Introduction. 


Many  problems  in  operations  research  involve  minim  ieing  a  cost  function  defined  on  a  bipartite 
or  directed  graph.  A  simple  but  fundamental  example  is  the  assignment  problem.  This  paper 
gives  algorithms  for  such  problems  that  run  almost  as  fast  as  the  best-known  algorithms  for  the 
corresponding  problems  without  costs.  For  the  assignment  problem,  the  corresponding  problem 
without  costs  is  maximum  cardinality  bipartite  matching. 

The  results  are  achieved  by  scaling  the  costs.  This  requires  the  costs  to  be  integral-valued. 
Farther,  for  the  algorithms  to  be  efficient,  costs  should  be  polynomially-bounded  in  the  number 
of  vertices,  >.e.,  n°l^.  These  requirements  are  satisfied  by  a  large  number  of  problems  in  both 
theoretical  and  practical  applications. 

Table  I  summarizes  the  results  of  the  paper.  The  parameters  describing  the  input  are  specified 
in  the  caption  and  defined  more  precisely  below.  The  first  column  gives  the  problem  and  the 
best-known  strong-polynomial  time  bound.  Such  a  bound  comes  from  an  algorithm  with  run  time 
independent  of  the  size  of  the  numbers  (assuming  the  uniform  cost  model  of  computation  [AHU]). 
The  second  column  gives  the  time  bounds  achieved  in  this  paper  by  scaling.  The  table  shows 
that  significant  speed-ups  can  be  achieved  through  scaling.  Further  it  will  be  seen  that  the  scaling 
algorithms  are  simple  to  program.  Now  we  discuss  the  specific  results. 

The  assignment  problem  is  to  find  a  minimum  cost  perfect  matching  on  a  bipartite  graph.  The 
strong  polynomial  algorithm  is  the  Hungarian  algorithm  [K55,K56]  implemented  with  Fibonacci 
heaps  (FT].  This  algorithm  can  be  significantly  improved  when  all  costs  are  zero.  Then  the  problem 
amounts  to  finding  a  perfect  matching.  The  best-known  cardinality  matching  algorithm,  due  to 
Hopcroft  and  Karp,  runs  in  time  0(y/nm)  (HKj.  The  new  time  bound  for  the  assignment  problem 
is  just  a  factor  lognJV  more  than  this.  The  algorithm  is  similar  to  the  Hopcroft-Karp  cardinality 
algorithm  and  appears  simple  enough  to  be  useful  in  practice. 

The  new  algorithm  improves  the  scaling  algorithm  of  [G85a],  which  runs  in  time 
0(n9/4m  log  N).  The  improvement  comes  from  a  different  scaling  method.  The  algorithms  of  [G85a] 
compute  an  optimum  solution  at  each  of  log  N  scales.  The  new  method  computes  an  approximate 
optimum  at  each  of  lognA^  scales;  using  log  n  extra  scales  ensures  that  the  last  approximate  op¬ 
timum  is  exact.  The  appropriate  definition  of  approximate  optimum  is  due  to  Bertsekas  [Ber86]. 
The  new  approach  to  scaling  was  recently  discovered  by  Goldberg  and  Tarjan  for  the  minimum  cost 
flow  problem  (Go,  GoT87a-b].  Their  minimum  cost  flow  algorithm  solves  the  assignment  problem 
in  time  0(nm  log(niV)),  which  this  paper  improves.  Bertsekas  [Ber87]  gives  an  algorithm  for  the 


assignment  problem  that  also  runs  sequential  time  0(nm  log (n IV)),  and  in  practice  appears  to  run 
faster  on  parallel  machines. 

The  assignment  algorithm  extends  to  other  network  problems.  Variants  of  minimum  cost 
perfect  matching  (such  as  miniwinm  cost  matching)  can  be  done  in  the  same  time  bound.  The 
linear  programming  dual  variables  for  perfect  matching  can  be  found  from  the  algorithm.  This 
gives  a  solution  to  the  shortest  path  problem  when  negative  edge  lengths  are  allowed.  The  table 
entry  for  the  degree-constrained  subgraph  problem  is  just  a  factor  log  nN  more  than  the  bound  of 
[ET]  for  the  corresponding  problem  without  costs,  namely  the  problem  of  maximum  flow  on  a  0*1 
network.  These  bounds  improve  (G85aj  in  a  manner  analogous  to  matching. 

The  table  entry  for  the  transportation  problem  is  a  good  bound  when  total  supply  and  demand 
(U)  is  small.  The  key  fact  for  this  bound  is  the  low  total  augmenting  path  length  for  the  assignment 
algorithm;  this  fact  generalizes  the  bounds  of  (ET]  for  cardinality  matching  and  0*1  network  flow. 
The  entry  for  minimum  cost  flow  is  a  double  scaling  algorithm—  it  scales  edge  capacities,  and  at  each 
scale  solves  a  small  transportation  problem  by  the  above  cost-scaling  algorithm.  This  algorithm  is 
not  as  good  asymptotically  as  the  recent  bound  of  Goldberg  and  Tarjan,  0(nmlog(n2/m)  log(nJV)) 
[GoT87b].  The  latter  is  just  a  factor  log  nN  more  than  the  best  bound  for  maximum  value  flow 
(G0T86].  However  the  double  scaling  algorithm  may  be  more  useful  in  practice,  since  it  requires 
fewer  data  structures.  Furthermore  it  generalizes  to  allow  the  cost  of  an  edge  to  be  an  arbitrary 
convex  function  of  the  flow  in  the  edge;  the  time  bound  is  unchanged,  as  long  as  the  cost  for  a  given 
flow  value  can  be  computed  in  0(1)  time.  The  last  table  entry,  for  the  Chinese  postman  problem, 
is  a  consequence  of  the  mincost  flow  algorithm.  The  algorithm  scales  capacities  but  the  bound  is 
strongly  polynomial. 

Section  2  presents  the  matching  algorithm  and  its  analysis,  including  facts  used  in  the  general¬ 
izations.  Section  3  presents  the  extensions  to  more  general  network  flow  problems.  Section  4  gives 
brief  concluding  remarks.  This  section  closes  with  definitions  from  graph  theory;  more  thorough 
treatments  are  in  (L,  T83]. 

We  use  interval  notation  for  sets  of  integers.  Thus  for  integers  i  and  j, 

[»../]  *  {k\k  is  an  integer, »  <  j  <  k),  (»../)  «  {k\k  is  an  integer, »  <  j  <  Jb},  etc.  The  function  log 
n  denotes  logarithm  to  the  base  two. 

For  a  graph  G,  V(G)  and  E(G)  denote  the  vertex  set  and  edge  set,  respectively.  The  given 
graph  G  is  bipartite  and  has  bipartition  V0,  V\  (so  V(G)  partitions  into  Vo  and  Vi,  and  any  edge 
joins  Vo  to  Vj).  The  given  graph  G  has  m  edges;  in  Section  2,  n  =  |Vo|  *  |Vj|  (we  assume 
without  loss  of  generality  that  the  two  sets  of  the  bipartition  have  equal  cardinality);  in  Section  3, 


n  s  |V(G)|.  If  E  is  a  subgraph  of  G,  an  H-edge  is  an  edge  in  H  and  a  non-H-edge  is  not  in  H.  When 
an  auxilary  graph  G'  is  constructed  from  the  given  graph  G,  G-edge  refers  to  an  edge  of  G'  that 
represents  an  edge  of  G.  We  use  this  term  without  explicit  comment  only  when  the  representation 
is  obvious  (vw  €  E(G)  is  represented  by  vV  where  v1  and  w'  are  obvious  representatives  of  v  and 
w).  We  say  path  P  ends  with  edge  w  if  vw  is  at  an  end  of  P  and  further,  v  is  an  endpoint  of  P. 

A  matching  on  a  graph  is  a  set  of  vertex-disjoint  edges.  Thus  a  vertex  v  is  in  at  most  one 
matched  edge  vv1;  v1  is  the  mate  of  v.  A  free  vertex  has  no  mate.  A  maximum  cardinality  matching 
has  the  greatest  number  of  edges  possible;  a  perfect  matching  has  no  free  vertices  (and  is  clearly 
maximum  cardinality).  An  alternating  path  (cycle)  for  a  matching  is  a  simple  path  (cycle)  whose 
edges  are  alternately  matched  and  unmatched.  An  augmenting  path  P  is  an  alternating  path  joining 
two  distinct  free  vertices.  Augmenting  the  matching  along  P  means  enlarging  the  matching  M  to 
A/©P,  thus  giving  a  matching  with  one  more  edge.  Suppose  each  edge  e  has  a  numeric  cost  c(e);  in 
this  paper  costs  are  integers  in  [-  W..WJ  unless  6tated  otherwise.  The  cost  c(S)  of  a  set  of  edges  5  is 
the  sum  of  the  individual  edge  costs.  A  minimum  (maximum)  perfect  matching  is  a  perfect  matching 
of  smallest  (largest)  possible  cost.  The  assignment  problem  is  to  find  a  minimum  perfect  matching 
on  a  bipartite  graph.  More  generally,  a  minimum  cost  maztmum  cardinality  matching  is  a  matching 
that  has  the  greatest  number  of  edges  possible,  and  subject  to  that  restriction  has  minimum  cost 
possible.  (The  phrase  ‘minimum  cost  maximum  cardinality  set*  can  be  interpreted  ambiguously. 
In  this  paper  it  refers  to  a  set  that  has  maximum  cardinality  subject  to  any  other  restrictions  that 
have  been  mentioned,  and  among  such  sets  has  minimum  cost  possible).  A  minimum  cost  matching 
is  a  matching  of  minimum  cost  (its  cardinality  can  be  any  value,  including  zero). 

A  multigraph  has  a  set  of  edges  £(G),  where  each  edge  e  has  an  integral  multiplicity  u(e)  (i.e., 
there  are  u(e)  parallel  copies  of  c).  The  size  parameter  m  is  the  number  of  edges,  m  =  |£(G)|;  m 
counts  multiplicities,  i.e.,  171  =  £{u(e)|e  €  £((?)};  AT  is  the  maximum  edge  multiplicity.  (In  a  graph 
Ad  -  1).  When  each  vertex  v  has  associated  nonnegative  integers  £(v)  and  «(v),  a  degree- constrained 
subgraph  ( DCS)  is  a  subgraph  where  each  vertex  has  degree  in  [£(v)..u(v)J.  It  is  convenient  to  use 
both  set  notation  and  functional  notation  for  a  DCS.  Thus  we  use  a  capital  letter  D  to  denote  a 
DCS,  and  the  corresponding  lower  case  letter  d  to  denote  two  functions  defined  by  D:  for  an  edge 
e,  d{e)  denotes  the  multiplicity  of  e  in  D,  and  for  a  vertex  v,  d(v)  denotes  the  degree  of  v  in  D,  i.e., 
d(v)  *  ][){d(vw)|vt0  €  £(C)}.  Hence  d{e)  <  u(e)  and  t(v)  <  d(e)  <  u(v).  The  deficiency  of  DCS 
D  at  vertex  v  is  ^(v,  D)  =  v(v)  -  d(v).  In  a  perfect  DCS  each  deficiency  is  zero.  The  size  of  the 
DCS  is  measured  by  U  *  £{u(v)|v  €  V)  (so  V  is  twice  the  number  of  edges  in  a  perfect  DCS). 
When  edges  e  have  costs  the  usual  assumption  is  that  each  copy  of  t  has  the  same  cost,  denoted 
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c(e).  When  this  assumption  fails  we  use  cost  functions,  defined  in  the  text.  Other  definitions  for 
DCS —  e.g minimum  perfect  DCS ,  minimum  coat  maximum  cardinality  DCS,  minimum  cost  DCS, 
follow  by  analogy  with  matching. 

The  transportation  problem  is  to  find  a  mmimnTn  cost  perfect  DCS  in  a  bipartite  multigraph 
where  all  edges  have  infinite  multiplicity;  alternatively  if  M  is  the  maximum  degTee  constraint,  all 
multiplicities  are  M.  If  some  multiplicities  are  less  than  AS  the  problem  is  a  capacitated  transporta¬ 
tion  problem.  (The  usual  definition  of  the  transportation  problem  allows  nonnegative  real-valued 
degree  constraints  and  edge  multiplicities  (both  given  multiplicities  and  those  in  the  solution).  This 
paper  deals  with  the  integral  case  of  this  problem.  Note  that  if  the  given  degree  constraints  and 
multiplicities  are  rational,  they  can  be  scaled  up  to  integers.  Also  note  that  no  loss  of  generality 
results  from  the  constraint  in  this  paper  that  the  solution  to  the  transportation  problem  has  in¬ 
tegral  multiplicities —  such  an  optimum  solution  always  exists  when  the  given  degree  constraints 
and  multiplicities  are  integral  (L).  Finally  note  that  in  our  terminology  the  minimum  perfect  DCS 
problem  is  the  same  as  the  capacitated  transportation  problem). 

Next  consider  a  directed  graph  G  with  nonnegative  edge  costs.  The  directed  Chinese  post¬ 
man  problem  is  to  find  a  multiplicity  for  each  edge  of  G  so  the  corresponding  multigraph  M  is 
Eulerian,  or  if  this  is  not  possible  M  decomposes  into  the  fewest  possible  number  of  open  paths, 
and  subject  to  this  restriction,  M  has  the  smallest  cost  possible  [PS].  This  amounts  to  the  fol¬ 
lowing  problem:  Let  indegree(v)  and  ovtdegrcc{v)  be  the  number  of  edges  directed  to  v  and 
from  v,  respectively.  Call  vertex  v  a  start  vertex  if  indegree(v)  >  outdegree(v),  an  end  vertex 
if  indegree(v)  <  outdegree(v),  and  otherwise  an  intermediate  vertex.  The  problem  is  to  find  a 
minimum  cost,  maximum  cardinality  set  of  paths  P  that  start  only  at  start  vertices,  end  only  at 
end  vertices,  each  start  vertex  v  starts  at  most  indegree(v)  -  outdegree(v)  paths,  and  each  end 
vertex  v  ends  at  most  outdegree(v)  -  indcgree(v)  paths.  (The  given  graph  G  becomes  an  Eule¬ 


rian  multigraph  if  P  contains  p  paths,  for  p  =  ^{indegree(v)  -  ovtdegree(v) |v  a  start  vertex)  = 
£{outdegree( v)  -  *Weyree(v)|v  an  end  vertex}). 

2.  Matching  and  extensions. 

Section  2.1  presents  our  algorithm  to  find  a  minimum  perfect  matching  on  a  bipartite  graph. 
Section  2.2  gives  extensions  to  other  versions  of  matching,  some  facts  about  the  algorithm  needed 
in  Section  3,  and  our  shortest  path  algorithm.  In  this  section  n  denotes  the  number  of  vertices  in 
each  vertex  set  V0,  Vt  of  the  given  bipartite  graph. 
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2.1.  The  matching  algorithm. 

For  convenience  assume  the  given  graph  G  has  a  perfect  matching  (the  algorithm  can  detect 
graphs  not  having  a  perfect  matching,  as  indicated  below). 

The  plan  for  the  algorithm  is  to  combine  the  Hungarian  algorithm  for  weighted  matching  with 
the  Hopcroft-Karp  algorithm  for  cardinality  matching.  Recall  that  the  Hungarian  algorithm  always 
chooses  an  augmenting  path  of  smallest  net  cost.  The  Hopcroft-Karp  algorithm  always  chooses  an 
augmenting  path  of  shortest  length.  Both  of  these  rules  can  be  approximated  simultaneously  if 
the  costs  are  small  integers.  Arbitrary  costs  can  be  replaced  by  small  integers  by  scaling.  Thus 
our  algorithm  scales  the  costs.  At  each  scale  it  computes  a  matching.  The  computation  is  efficient 
because  it  is  similar  to  the  Hopcroft-Karp  algorithm;  the  matching  is  close  to  optimum  because  the 
computation  is  similar  to  the  Hungarian  algorithm.  Now  we  give  the  details. 

Each  scale  of  the  algorithm  finds  a  close-to-minimum  matching,  defined  as  follows.  Every 
vertex  v  has  a  dual  variable  y(v).  A  1 -feasible  matching  consists  of  a  matching  M  and  dual 
variables  y(t>)  such  that  for  any  edge  vw, 

y(«)  4-  y(u>)  <  e(i>to)  4-  1, 

y(t>)  4-  y(w)  *  c(vw),  for  vid  6  M. 

A  1~  optimal  matching  is  a  perfect  matching  that  is  1-feasible.  If  the  4-1  term  is  omitted  horn 
the  first  inequality,  these  are  the  usual  complementary  slackness  conditions  for  a  minimum  perfect 
matching  (L).  The  following  result  is  due  to  Bertsekas  [BerS6j. 

Lemma  2.1.  Let  M  be  a  1 -optima)  matching. 

(а)  Any  perfect  matching  P  has  c(P)  >  c[M)  -  n. 

(б)  If  some  integer  k,  k  >  n,  divides  each  cost  c(e)  then  M  is  a  minimum  perfect  matching. 

Proof.  Part  (a)  follows  because  c[M)  -  £{c(e)|e  €  M }  =  £{y(v)|v  €  V^G)}  <  c(P)  4-n.  Part 
(6)  follows  from  (a)  and  the  fact  that  any  matching  has  cost  a  multiple  of  k.  I 

This  lemma  is  the  basis  for  the  main  routine  of  the  algorithm,  which  does  the  scaling.  The 
routine  starts  by  computing  a  new  cost  ?(e)  for  each  edge  e,  equal  to  n  4- 1  times  the  given  cost. 
Consider  each  7(e)  to  be  a  signed  binary  number  ±6162  ■  •  ■  bk  of  k  =  [log(n  4-  1)1VJ  4- 1  bits.  The 
routine  maintains  a  variable  c(e)  for  each  edge  e,  equal  to  its  cost  in  the  current  scale.  The  routine 
initializes  each  c(e)  to  0  and  each  dual  y(v)  to  0.  Then  it  executes  the  following  loop  for  index  t 
going  from  1  to  k: 
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Step  1.  For  each  edge  e,  e(e)  •-  2c(e)+  (signed  bit  b,  of  7(e)).  For  each  vertex  v,  y(v)  —  2y(v)  - 1. 
Step  2.  Call  the  scale jmatck  routine  to  find  a  1  •optimal  matching.  I 

Lemma  2.1(b)  shows  that  the  routine  halts  with  a  minimum  perfect  matching.  Each  iteration 
of  the  loop  is  called  a  scale.  We  give  a  scalejmatch  routine  that  runs  b  O(^nm)  time.  Sbce  there 
are  0(log(nAr))  scales  this  achieves  the  desired  time  bound. 

It  is  most  natural  to  work  with  small  costs.  The  scalejmatch  routine  transforms  costs  to  achieve 
this.  Specifically,  scale-match  changes  the  cost  of  each  edge  vw  to  c(vw)  -  y(v)  -  y(w)\  then  it  calls 
the  match  routbe  on  these  costs  to  find  a  1 -optimal  matchbg  M  with  duals  y'(t');  then  it  adds 
y*(v)  to  each  dual  y(v)  (y(v)  is  the  dual  value  before  the  call  to  match). 

Clearly  M  with  the  new  duals  is  a  1 -optimal  matchbg  for  cost  function  c.  Further,  sbce  Step 
1  of  the  mab  routbe  changes  costs  and  duals  so  that  the  empty  matchbg  is  1 -feasible,  the  costs 
bput  to  match  are  btegers  -1  or  larger.  If  vw  is  an  edge  b  the  1 -optimal  matchbg  found  b  the 
previous  scale  then  after  Step  1,  y[v)  +  y(w)  >  c(vw)-  3.  Hence  vw  costs  at  most  three  b  the  costs 
for  match.  Thus  there  is  a  matchbg  of  C06t  at  most  3n.  (This  is  true  even  b  the  first  scale).  We 
will  show  that  if  every  edge  costs  at  least  -1  and  a  minimum  perfect  matching  costs  O(n),  match 
finds  a  1 -optimal  matchbg  b  0(y/nm)  time.  This  gives  the  desired  time  bound. 

Note  that  the  transformation  done  by  scalejmatch  is  for  conceptual  convenience  only.  An  actual 
implementation  would  not  transform  costs;  rather  match  would  work  directly  on  the  untransformed 
costs. 

The  cost-length  of  an  edge  e  with  respect  to  a  matchbg  M  is 

ef(e)  =  c(e)  +  (if  e  M  then  1  else  0). 

The  net  cost-length  of  a  set  of  edges  5  with  respect  to  M  is 

e/(S)  *  £{c/(e)|e  £  S  -  M)  -  £{e/(e)je  €  S  n  M). 

This  quantity  equals  the  net  cost  of  5  (with  respect  to  M)  plus  the  number  of  unmatched  edges 
b  S.  Hence  an  augmentbg  path  with  smallest  net  cost-length  approximates  both  the  smallest  net 
cost  augmentbg  path  and  the  shortest  augmentbg  path;  this  is  b  keeping  with  our  plan  for  the 
algorithm. 

An  edge  vw  is  eligible  if  y(v)+y(w)  =  cl(vw),  i.e.,  the  1-feasibility  constraint  for  vw  holds  with 
equality.  (Note  that  a  matched  edge  is  always  eligible).  It  follows  from  the  analysis  below  that  an 
augmentbg  path  of  eligible  edges  has  the  smallest  possible  net  cost-length.  Hence  the  algorithm 


augments  along  paths  of  eligible  edges.  If  no  such  path  exists  it  adjusts  the  duals  to  create  one. 
The  details  are  as  follows. 

Assume  the  costs  given  to  match  are  integers  at  least  —1,  and  there  is  a  perfect  matching 
costing  at  most  an.  (In  the  scaling  algorithm  a  =  3). 

procedure  match. 

Initialize  all  duals  y(v)  to  0  and  matching  A/  to  0.  Then  repeat  the  following  steps  until  Step  1 
halts  with  the  desired  matching: 

Step  1.  Find  a  maximal  set  A  of  vertex  disjoint  augmenting  paths  of  eligible  edges.  For  each  path 
P  €  A,  augment  the  matching  along  P,  and  for  each  vertex  w  €  V\  n  P,  decrease  y(w)  by  1.  (This 
makes  the  new  matching  A/  1 -feasible).  If  the  new  matching  U  is  perfect,  halt. 

Step  2.  Do  a  Hungarian  search  to  adjust  the  duals  (maintaining  1-feasibility)  and  find  an  aug¬ 
menting  path  of  eligible  edges.  I 

We  now  give  the  details  of  Steps  1  and  2  that  are  needed  to  analyze  match.  (A  full  description 
of  these  steps  is  given  later).  Both  steps  can  be  implemented  in  O(m)  time.  Step  2  is  a  Hungarian 
search  (essentially  Dijkstra’s  shortest  path  algorithm,  see  e.g.  [L]).  The  search  does  a  number  of 
dual  adjuttmenti.  Each  dual  adjustment  calculates  a  positive  integer  6  and  increases  or  decreases 
various  dual  values  by  6,  so  as  to  preserve  1-feasibility  and  eventually  create  an  augmenting  path 
of  eligible  edges.  (The  dual  adjustment  is  defined  more  precisely  below).  At  any  point  in  the 
algorithm  define 

F  =  the  set  of  free  vertices  in  Vo; 

A  —  the  sum  of  all  dual  adjustment  quantities  6 
in  all  Hungarian  searches  so  far. 

The  Hungarian  search  maintains  the  duals  so  that  any  free  vertex  v  €  F  has  y(v)  -  A  and  any 
free  vertex  v  €  V'i  has  y(v)  =  0. 

To  analyze  the  match  routine,  first  observe  that  it  is  correct:  The  changes  to  the  matching  (in 
Step  l)  and  to  the  duals  (in  Step  2)  keep  A/  a  1 -feasible  matching.  If  A/  is  not  perfect  but  G  has 
a  perfect  matching,  the  Hungarian  search  creates  an  augmenting  path  of  eligible  edges.  Hence  the 
algorithm  eventually  halts  with  M  a  1-optimal  matching,  as  desired.  (Note  that  if  G  does  not  have 
a  perfect  matching,  this  is  eventually  detected  in  Step  2). 

To  analyze  the  run  time,  consider  any  point  in  the  execution  of  match.  Let  M  be  the  current 
matching,  and  define  F  and  A  as  above.  Let  M*  be  a  minimum  perfect  matching.  For  any  set  of 
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To  analyze  the  run  time,  consider  any  point  in  the  execution  of  match.  Let  M  be  the  current 
matching,  and  define  F  and  A  as  above.  Let  M*  be  a  minimum  perfect  matching.  For  any  set  of 
edges  S  let  cl(S)  denote  net  cost-length  with  respect  to  M. 

M‘  ©  M  consists  of  an  augmenting  path  P,  for  each  tiff,  plus  alternating  cycles  Cw.  Thus 

n  +  c(AT)  -  c(M)  >  c/(Af  ©  M)  =  g c/(P.)  +  £  c/(C„).  (1) 

To  estimate  the  right-hand  side,  consider  an  alternating  path  P  from  u  €  Vo  to  m  €  Vo,  where 
tc  is  on  an  unmatched  edge  of  P  and  m  is  on  a  matched  edge  of  P  (m  stands  for  ‘matched”;  no 
confusion  should  result  from  the  double  usage  of  m).  Then 

y(u)  <  y(m)  +  c/(P).  (2) 

This  follows  since  for  edges  vv$M  and  vm  €  M,  y(«)+y(v)  <  cl(uv)  and  y(v)+y(m)  =  cl(vra),  so 
y(u)  <  y[m)  +  c/(uv)  —  cl(vm).  Inequality  (2)  implies  that  any  alternating  cycle  C  has  cl[C)  >  0. 
It  also  implies  that  any  augmenting  path  P„  from  some  v  €  F  to  some  free  vertex  t  €  V\  has 
y(r)  +  y(t)  <  c/(Pv).  Recall  that  the  Hungarian  search  keeps  y(v)  =  A  and  y(t)  =  0.  Hence 
A  <  c/(Pv),  and  the  right-hand  side  of  (1)  is  at  least  |F|A. 

By  assumption  on  the  input  to  match,  c(M”)  <  an  and  c(M)  >  -n.  Hence  the  left-hand  side 
of  (1)  is  at  most  6n  for  b  «  a  +  2.  Thus  we  have  shown 

|F|A  <  bn.  (3) 

This  implies  there  are  0(y/n)  iterations  of  the  loop  of  match.  To  see  this  note  that  each 
execution  of  Step  1  (except  possibly  the  first)  augments  along  at  least  one  path,  because  of  the 
preceding  Hungarian  search.  Hence  at  most  y/Fn  + 1  iterations  start  with  |F|  <  y/Fn.  From  (3), 
|F|  >  y/Fn  implies  A  <  y/Fn.  We  will  show  that  each  Hungarian  search  increases  A  by  at  least 
one.  This  implies  at  most  y/Fn  +  1  iterations  start  with  A  <  y/Fn,  giving  the  desired  bound. 

Now  we  show  that  a  Hungarian  search  S  increases  A  by  at  least  one.  It  suffices  to  show  that  5 
does  a  dual  adjustment  (since  any  dual  adjustment  quantity  6  is  a  positive  integer).  Search  S  does 
a  dual  adjustment  unless  when  it  starts,  there  is  an  augmenting  path  P  of  eligible  edges.  Clearly  P 
intersects  some  augmenting  path  found  in  Step  1.  It  is  easy  to  see  that  P  contains  an  unmatched 
edge  vw,  with  w  but  not  v  in  an  augmenting  path  of  Step  1,  and  w  6  Vj .  But  when  Step  1  decreases 

y(u>)  this  makes  vw  ineligible.  So  P  does  not  exist,  and  5  does  a  dual  adjustment. 

In  summary  match  does  0(v/n)  iterations.  Each  iteration  takes  time  O(m),  giving  the  desired 
time  bound  0{y/nm). 
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It  remains  to  give  the  details  of  Steps  1  and  2.  Step  1  finds  the  augmenting  paths  P  by  depth- 
first  search.  To  do  this  it  marks  every  vertex  reached  in  the  search.  It  initialises  a  path  P  to  a  free 
unmarked  vertex  of  Vo.  To  grow  P  it  scans  an  eligible  edge  xy  from  the  last  vertex  x  oi  P  (x  will 
always  be  in  Vo).  If  y  is  marked,  the  next  eligible  edge  from  z  is  scanned;  if  none  exists  the  last  two 
edges  of  P  (one  matched  and  one  unmatched)  are  deleted  from  P\  if  P  has  no  edges  another  path 
is  initialized.  If  y  is  free  another  augmenting  path  has  been  found,  y  is  marked,  the  path  is  added 
to  A  and  the  next  path  is  initialized.  The  remaining  possibility  is  that  y  is  matched  to  a  vertex  z. 
In  this  case  y  and  z  are  marked,  edges  zy,  yz  are  added  to  P  and  the  search  is  contbued  from  z. 

It  is  clear  that  this  search  uses  0(m)  time.  To  show  that  it  halts  with  A  maximal,  first  observe 
that  for  any  marked  vertex  z  €  V0  -  V(A),  every  eligible  edge  xy  has  y  marked  and  matched,  or  b 
V(A).  (Note  that  V(A)  is  the  set  of  vertices  b  paths  of  A).  Hence  an  easy  bduction  shows  that 
an  alternatbg  path  of  eligible  edges,  startbg  at  a  free  vertex  of  Vo  and  vertex  disjobt  from  A ,  has 
all  its  Vo  vertices  marked,  and  is  not  augmentbg. 

Step  2  is  the  Hungarian  search.  It  grows  a  forest  7  of  eligible  edges,  from  roots  F.  A  pair  of 
eligible  edges  vte,  ww',  where  v  €  V0  n  7,  w  ^  7,  and  ww'  6  M,  is  added  to  7  whenever  possible. 
Eventually  either  an  augmentbg  path  of  eligible  edges  is  found  or  7  cannot  be  enlarged. 

In  the  latter  case  a  dual  adjustment  is  done.  It  changes  duals  b  a  way  that  preserves  1- 
feasibility  and  allows  7  to  be  enlarged,  as  follows.  It  computes  the  dual  adjustment  quantity 

6  s=  mb{c/(vty)  -  y(v)  -  y(u>)|v  €  Vo  n  7,  7). 

For  each  v  6  7,  it  bcreases  y(v)  by  6x  (if  v  €  Vo  then  1  else  —1).  It  is  easy  to  see  that  this 
achieves  the  goal  of  the  dual  adjustment  (any  edge  vw  achievbg  the  above  minimum  becomes 
eligible  after  the  dual  adjustment,  and  so  can  be  added  to  7). 

After  the  dual  adjustment  the  search  contbues  by  enlargbg  7.  Eventually  the  desired  aug¬ 
mentbg  path  is  found. 

Note  that  as  claimed  above,  at  any  pobt  b  the  algorithm  a  free  vertex  v  has  y(v)  =  A  if  v  €  F 
(sbce  every  dual  adjustment  bcreases  y(v))  and  y(v)  =  0  if  v  €  V\  (no  dnal  adjustment  changes 

*(»)) 

A  Hungarian  search  can  be  implemented  b  0(m)  time.  This  depends  on  two  observations. 
First,  the  proper  data  structure  allows  a  dual  adjustment  to  change  all  duals  y(v)  in  0(1)  time 
total.  Specifically  the  algorithm  keeps  track  of  A  (defined  above).  When  a  vertex  t>  is  added  to  7 
its  current  dual  value  and  the  current  value  of  A  are  saved  as  y°(v)  and  A(v),  respectively.  Then 
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at  any  time  the  current  value  of  y(v)  can  be  calculated  as 

y°(t>)  +  (A  —  A(t>))  x  (If  v  €  Vo  then  1  else  -  1). 

Hence  the  dual  adjustment  is  accomplished  by  simply  increasing  the  value  of  A. 

The  second  observation  is  how  to  compute  £  in  a  dual  adjustment.  The  usual  implementation 
of  a  Hungarian  search  does  this  with  a  priority  queue  that  introduces  a  logarithmic  factor  into  the 
time  bound  [e.g.,  FT].  This  can  be  avoided  when,  as  in  our  case,  costs  are  small  integers  (this  was 
observed  in  [D,  W]  for  Djjkstra’s  shortest  path  algorithm).  The  details  are  as  follows.  The  next 
value  of  6  is  the  amount  that  the  next  value  of  A  increases  from  its  current  value.  Hence  it  suffices 
to  calculate  the  next  value  of  A.  The  next  value  of  A  is  the  smallest  possible  value  such  that  some 
edge  vw  with  v  €  Vo  fl  7  and  w  7  becomes  eligible  (when  duals  are  adjusted  by  6).  Thus  the 
next  value  of  A  equals 

min{e!(vw)  -  y°(t>)  -  y(w)  +  A(ti)|v  €  Vo  n  7,  w^J). 

Since  any  Hungarian  search  has  |F|  >  1,  inequality  (3)  implies  A  <  bn.  The  algorithm  maintains 
an  array  <?[l..bn].  Each  entry  Q[r]  points  to  a  list  of  edges  vw  that  can  make  A  =  r,  i.e.,  v  € 
V0  n  7,  w  t  7  and  r  *  el[vw)  -  y°(v)  -  y{w)  +  A(v).  The  algorithm  scans  down  Q  and  chooses 
next  value  of  A  as  the  smallest  value  r  with  Q\r]  nonempty.  This  gives  the  next  value  of  6,  and  the 
newly  eligible  edges,  as  desired.  The  total  overhead  for  scanning  is  O(n)  since  Q  has  bn  entries. 
(Note  that  an  edge  vw  with  v  6  V0  n  7,  w  7  does  not  get  entered  in  this  data  structure  if 
el(vw)  -  y°(«)  -  y(u>)  +  A(v)  >  bn). 

Only  one  detail  of  the  derivation  remains:  We  have  assumed  that  the  dual  values  y(v)  do  not 
grow  too  large,  so  that  arithmetic  operations  use  0(1)  time.  To  justify  this  we  show  that  each  y(v) 
has  magnitude  0(n2N).  It  suffices  to  do  this  for  v  €  Vo.  Define  Y,  as  max(|y(v)||v  €  V0)  after  the 
8tk  scale.  Then  Vo  =  Oand  VVh  <  2V,+bn-l  (since  A  <  bn).  ThusV*  <  (2*-l)(bn-l)  =  0(n2N) 
as  desired.  Note  that  the  input  uses  a  word  sue  of  at  least  max{iog  N,  logn)  bits.  Bence  at  worst 
the  algorithm  uses  triple-word  integers  for  the  dual  variables. 

Theorem  2.1.  A  minimum  perfect  matching  on  a  bipartite  graph  can  be  found  in 
0(>/nmlog(nAf))  time  and  0(m)  space.  I 

A  heuristic  that  may  speed  up  the  algorithm  in  practice  is  to  prune  the  graph  at  the  start  of 
each  scale.  Specifically,  scale-match  can  delete  any  edge  whose  new  cost  is  6n  or  more.  In  proof, 
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recall  that  in  the  costs  computed  by  acalcjmatch  there  is  a  perfect  matching  M  costing  at  most  3n; 
taking  into  account  the  low  order  bits  of  cost  that  are  not  included  in  the  current  scale,  the  true 
cost  of  M  is  less  than  in.  In  the  costs  computed  by  scale. match  every  edge  costs  at  least  -1;  again 
taking  into  account  low  order  bits,  the  true  cost  is  more  than  —2.  Hence  a  matching  containing  an 
edge  of  new  cost  6n  or  more  has  true  cost  more  than  4n  and  so  is  not  minimum. 

2.2.  Extensions  of  the  matching  algorithm. 

The  bounds  of  Theorem  2.1  also  apply  to  finding  a  minimum  cost  matching.  To  see  this  let  G 
be  the  given  graph.  Form  V  by  taking  two  copies  of  G;  for  each  v  €  V(G)  join  the  two  copies  of  v 
by  a  cost  zero  edge.  Then  G  is  bipartite,  and  a  minimum  perfect  matching  on  (7  gives  a  minimum 
cost  matching  on  G. 

A  similar  result  holds  for  minimum  cost  maximum  cardinality  matching.  The  construction  is 
the  same  except  that  the  edges  joining  two  copies  of  v  cost  nN .  The  problem  of  finding  a  minimum 
cost  matching  of  given  cardinality  can  also  be  solved  in  the  same  bounds;  it  is  most  convenient  to 
use  Theorem  3.2  below. 

Returning  to  perfect  matching,  several  properties  of  match  are  needed  for  Section  3.  Define 

A  =  the  total  length  of  all  augmenting  paths  found  by  match. 

We  first  derive  a  bound  on  A.  Let  Pi  denote  the  itk  augmenting  path  found  by  match.  Let  f, 
be  its  length,  measured  as  its  number  of  unmatched  edges;  let  A,  denote  the  value  of  A  when  jP, 
is  found;  let  Mx  be  the  matching  after  augmenting  along  Px.  Recall  that  in  the  Hopcroft-Karp 
algorithm,  for  some  constant  c,  f,  <  en/(n  —  i  +  1).  Thus  the  total  augmenting  path  length  is 
£"=1 —  0(n  log  n)  (ETJ.  In  match,  ti  does  not  have  a  similar  bound.  However  it  is  bounded  in 
an  amortized  sense,  as  follows. 

Lemma  2.2.  For  any  k  in  (l..n),  4  +  c(Afa)  =  £*ai  A,. 

Proof.  A  calculation  similar  to  (2)  shows  that  for  any  t,  A,  =  cl(Pi).  It  is  easy  to  see  that 
cl(Pi)  =  £,  +  c(A {{)  -  c(Af,-i)  (assume  c(Afb)  *  0).  Summing  these  relations  gives  the  lemma.  I 


Corollary  2.1.  A  =  0{n  log  n). 


Proof.  Since  |M*|  -  k  the  entry  conditions  for  match  imply  c(M*)  >  -k.  Hence  A  <  n  + 
£*Lj  A,.  By  (3),  A,  <  6n/(n  -  i  +  1).  Summing  these  inequalities  gives  the  lemma.  I 


The  second  property  shows  that  the  depth-first  search  of  Step  1  never  encounters  a  cycle.  A 
similar  property  for  network  flows  is  used  in  (GoT87a]. 


Lemma  2.3.  In  match  there  is  never  an  alternating  cycle  of  eligible  edges. 


Proof.  Initially  there  are  no  matched  edges,  so  there  are  no  alternating  cycles  of  eligible  edges. 
In  a  Hungarian  search,  whenever  the  duals  of  a  matched  edge  vw  are  changed,  w  6  Vj  gets  y(w) 
decreased.  Hence  any  edge  joining  w  to  a  vertex  not  in  the  search  forest  7  is  ineligible.  This  implies 
the  Hungarian  search  does  not  create  an  alternating  cycle  of  eligible  edges.  Similar  reasoning  applies 
when  an  augment  creates  a  new  matched  edge  and  changes  duals.  I 


Some  applications  of  matching  require  ordinary  dual  variables,  defined  as  follows.  The  duals 
are  dominated  on  edge  vw  if  y(v)  +  y(w)  <  c(vit>);  they  are  tight  on  vw  if  equality  holds.  The  duals 
are  dominated  and  tight  for  a  given  matching  if  each  edge  is  dominated  and  each  matched  edge  is 
tight;  such  duals  are  the  usual  linear  programming  dual  variables  [LJ.  The  scaling  algorithm  halts 
with  l-optimal  duals,  but  these  are  not  necessarily  dominated  and  tight.  Such  duals  can  be  found 
as  follows. 

Let  G+  be  G  with  an  additional  vertex  s  eV0  and  an  edge  »v  for  each  v  €  Vj .  Extend  the 
given  cost  function  c  to  G+  by  defining  c(sv)  as  an  arbitrary  integer;  the  cost  function  used  by  the 
matching  algorithm  extends  to  G+  by  its  definition,  ?  *  (n  +  l)c.  To  specify  a  cost  function  on  G+ 
we  write  <7+;c  or  G+;7.  Let  M  be  a  minimum  perfect  matching  on  G;  for  vertex  v  let  v'  denote 
its  mate,  Le.,  w'  €  M.  For  v  €  Vo  let  AC  be  a  minimum  perfect  matching  on  G+  -  v;c.  (Such  a 
matching  exists,  for  instance  M  —  vv'  +  at/).  Set 


y(v)  *  If  v  €  V0  then  -  c(Mv)  else  c(w')  -  y(v'). 


These  duals  are  dominated  and  tight  on  G.  (This  can  be  proved  by  an  argument  similar  to  the 
algorithm  given  below.  Alternatively  see  [G87]  for  a  proof  from  first  principles). 

Suppose  a  Hungarian  search  (as  in  match)  is  done  on  G+;7.  It  halts  with  a  tree  T  of  eligible 
edges,  rooted  at  s.  Clearly  T  is  a  spanning  tree.  For  any  v  6  V0,  augmenting  along  the  sv-path  in 
T  gives  a  l-optimal  matching  Nv  on  G+  -  v;7.  Nv  is  a  minimum  perfect  matching  on  G+  -  v;e. 
This  follows  from  Lemma  2.1,  since  G+  -  v  and  G  have  the  same  number  of  vertices.  Hence  N, 
qualifies  as  M,. 
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In  summary,  the  following  procedure  finds  dominated  tight  duals.  Given  is  the  output  of 
the  matching  algorithm,  i.e.,  a  1 -optimal  matching  on  G;7  with  duals  y.  Form  G+;r,  defining 
c(sv)  *  fy(v)/(n  -f  1)]  for  each  v  €  V\\  also  set  y(t)  *-  0  (this  gives  l-feasible  duals).  Do  a 
Hungarian  search  to  construct  a  spanning  tree  T  of  eligible  edges  rooted  at  s.  Do  a  depth-first 
search  of  T  to  find  e(. Mv)  for  each  v  6  Vo-  Define  dominated  tight  duals  y(v)  by  the  above  formula. 

The  time  for  this  algorithm  is  O(m).  This  is  clear  except  perhaps  for  the  time  for  the  Hungarian 
search.  The  choice  of  c(*t>)  ensures  that  A  <  n.  Hence,  as  in  match,  the  Hungarian  search  can  be 
implemented  using  an  array  Q.  This  gives  O(m)  time. 

Corollary  2.2.  Dominated  and  tight  duals  on  a  bipartite  graph  can  be  found  in 
0[y/nm\of(nN))  time  and  O(m)  space.  I 

This  implies  the  next  result.  Consider  a  directed  graph  with  n  vertices,  m  edges,  and  arbitrary 
(possibly  negative)  edge  lengths. 

Theorem  2.2.  The  single-source  shortest  path  problem  on  a  directed  graph  with  arbitrary  edge 
lengths  can  be  solved  in  Of^/nmlogfnAf))  time  and  O(m)  space. 

Proof.  This  problem  can  be  solved  by  finding  dominated  tight  duals  on  a  bipartite  graph  whose 
costs  are  the  edge  lengths,  and  then  running  Dijkstra’s  algorithm  [G85a].  I 

Obviously  the  same  bound  holds  for  0(y/n)  sources. 

3.  Degree-constrained  subgraphs  and  extensions. 

This  section  extends  the  matching  algorithm  to  derive  the  last  four  bounds  of  Table  I.  Section 
3.1  gives  an  algorithm  for  the  minimum  perfect  degree-constrained  subgraph  problem,  deriving 
time  bounds  for  finding  a  degree-constrained  subgraph  and  for  solving  the  transportation  problem. 
Section  3.2  discusses  scaling  edge  multiplicites,  which  improves  the  bounds  when  edge  multiplicities 
are  large.  Section  3.3  extends  the  results  to  network  flow.  Throughout  Section  3,  n  denotes  the 
number  of  vertices  in  the  input  graph. 

3.1.  The  degree-constrained  subgraph  algorithm. 
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This  subsection  gives  an  algorithm  for  the  perfect  degree-constrained  subgraph  problem.  Recall 
that  a  multigraph  has  m  edges  and  m  edges  counting  multiplicities.  Note  that  a  perfect  DCS 
problem  on  a  multigraph  of  n  vertices  and  Ifl  edges  can  be  reduced  in  linear  time  to  a  perfect 
matching  problem  on  a  graph  of  O(Tff)  vertices  and  edges  [G87].  Hence  Theorem  2.1  immediately 
implies  a  bound  of  Olflfit2  Jog(wrJV))  for  the  DCS  problem.  We  now  derive  the  better  bound  given 
in  Table  1. 

For  a  DCS  D,  the  coat-length  of  edge  e  is 

e/(e)  =  e(e)  +  (If  D  then  1  eke 0). 

A  1-feasible  DCS  is  a  DCS  D  and  dual  variables  y(v)  for  each  vertex  v,  such  that  for  any  edge  vw, 

y(v)  +  y(w)  <  cl(vw),  for  vw^D, 

y{v)  +  y(ty)  >  cl(vw),  for  vw  6  D. 

A  1-optimal  DCS  is  a  perfect  DCS  that  is  1 -feasible.  (Note  that  the  definition  of  a  1-feasible 
matching  is  slightly  different—  the  second  relation  holds  with  equality.  The  difference  is  not 
significant:  if  we  treat  a  matching  problem  as  a  DCS  problem,  a  1-feasible  DCS  gives  a  1 -feasible 
matching,  by  lowering  duals  as  necessary  to  achieve  the  desired  equalities). 

As  in  Lemma  2.1,  if  every  cost  is  divisible  by  k,  k  >  n/2,  then  a  1 -optimal  DCS  is  a  minimum 
perfect  DCS.  This  is  essentially  a  result  of  Bertsekas  [Ber86).  In  proof,  note  that  a  perfect  DCS  D 
has  minimum  cost  if  any  alternating  (simple)  cycle  C  has  c(C  D  D)  <  c{C  -  D).  This  condition 
can  be  verified  for  a  1-optimaf  DCS  D  by  a  calculation  similar  to  Lemma  2.1. 

Now  we  describe  the  algorithm.  Many  details  are  exactly  as  in  Section  2,  so  we  elaborate  only 
on  the  parts  that  change.  All  data  structures  have  size  0(m).  Clearly  the  multigraph  G  can  be 
represented  by  such  a  structure. 

The  main  routine  works  in  (at  most)  [log(n  +  2)AfJ  scales.  (This  is  justified  by  the  above 
analog  of  Lemma  2.1;  each  original  cost  is  multiplied  by  f(n  +  l)/2]).  Steps  1-2  and  scale-match 
are  unchanged.  Let  Do  be  the  1-optimal  DCS  of  the  previous  scale.  Note  that  the  match  routine 
is  called  with  integral  costs  c(e)  that  are  at  least  -1  for  e  £  Do  and  at  most  three  for  e  £  Do- 

The  match  routine  initializes  all  duals  y(v)  to  0  and  D  to  (e|c(e)  <  -1).  (Clearly  D  does  not 
violate  any  degree  constraint).  The  definition  of  an  eligible  edge  vw  is  still  y(v)  +  y(tc)  =  e/(vu>). 
Step  1  of  match  finds  a  maximal  set  of  augmenting  paths  of  eligible  edges  and  augments  the  matching 
along  each  path.  (Unlike  Section  2,  no  duals  are  changed  after  an  augment;  the  new  DCS  is  1- 
feasible  and  the  edges  on  an  augmenting  path  become  ineligible).  Step  2  does  a  Hungarian  search 
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to  adjust  duals  and  find  an  augmenting  path  of  eligible  edges.  Note  this  algorithm  is  cornet:  Since 
the  Hungarian  search  maintains  1-feasiblity,  the  algorithm  halts  with  a  1 -optimal  DCS  (assuming 
a  perfect  DCS  exists). 

Step  1  is  implemented  by  a  depth-first  search  similar  to  Section  2,  modified  for  degree  con¬ 
straints  larger  than  one:  Each  augmenting  path  P  is  initialised  to  a  vertex  z  £  VD  with  positive 
deficiency;  z  is  used  to  initialise  paths  P  until  its  deficiency  becomes  aero  or  it  is  deleted  from 
P.  Pa  grown  as  an  alternating  path,  so  that  when  its  last  vertex  z  is  in  Vo  an  edge  not  in  D  is 
scanned,  and  when  z  is  in  Vi  an  edge  of  D  is  scanned.  Instead  of  vertex  marks,  each  vertex  has 
a  pointer  to  its  last  unscanned  edge.  The  last  edge  of  P  gets  deleted  if  z  has  no  more  unscanned 
edges.  It  is  easy  to  see  the  time  for  Step  1  is  0(lfi).  (As  shown  below  each  augmenting  path  is 
simple,  although  this  fact  is  not  needed  for  correctness). 

The  details  of  the  Hungarian  search  are  similar  to  Section  2.  The  main  differences  stem  from 
the  fact  that  the  search  forest  7  is  grown  edge-by-edge,  rather  than  in  pairs  of  unmatched  and 
matched  edges.  The  time  for  the  search  is  O(m).  This  assumes  that,  as  in  Section  2,  an  array 
Q[l..dn]  is  used  to  compute  minima;  here  d  is  the  constant  of  Lemma  3.3,  which  justifies  using  the 
array. 

This  completes  the  description  of  the  DCS  algorithm.  The  discussion  shows  that  it  is  correct. 
The  efficiency  analysis  uses  three  inequalities,  each  analogous  to  (3)  of  Section  2.  We  use  notation 
similar  to  Section  2:  D  is  the  DCS  at  any  point  in  match.  Do  is  the  1-optimal  DCS  of  the  previous 
jeale;  hence  each  of  its  edges  costs  at  most  a  =  3.  F  is  the  set  of  vertices  in  V0  with  positive 
deficiency;  4>  is  their  total  deficiency, 

*-2>k  D)\veF). 

A  is  the  sum  of  all  dual  adjustment  quantities  £  in  all  Hungarian  searches  so  far.  Each  z  €  F  has 
y(z)  =  A.  P,  denotes  any  one  of  the  augmenting  paths  containing  z  in  Do  ©  D. 

Lemma  3.1.  For  some  constant  b,  at  any  point  in  match ,  $A  <  bV. 

Proof.  The  argument  of  Section  2  gives  an  analog  of  (1), 

cl(Do  ©  D)  >  D)y{v)\v  €  P). 

An  edge  of  Do  -  D  has  cost-length  at  most  a  •+ 1;  an  edge  of  D  -  Do  has  cost-length  at  least  -1. 
Hence  the  lemma  holds  with  6  =  (a  +  2)/2.  I 
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The  second  inequality  b  for  graphs  with  hounded  multiplicity.  It  generalises  [ETj.  Recall  that 
bd  denotes  the  maximum  multiplicity  of  an  edge  in  the  multigraph. 

Lemma  3.2.  For  some  constant  c,  at  any  point  in  match,  Av/$  <  cny/JH. 

Proof.  Set  6  =  a  4  2.  For  each  integer  j  define 

Ui~{ueVQ\y(u)z\b(j-\)..bj)l 

Wj  =  (w  €  Vijy(w)  -  s  -  1  €  (-6;..  -  6(j  -  1)]}. 

We  will  show  that  for  any  ;  €  (l..fA/6  4  1]],  each  P,  has  an  edge  uw  with  u  €  Uj,  w  6  Wr  This 
implies  M\Uj\\Wj\  >  4>.  Thus  |£/,|  or  \Wj\  is  at  least  y/¥fM.  Hence  n  >  y/$/M{&/b),  as  desired. 

To  find  the  desired  edge  uw  in  P, ,  let  the  edges  in  Pt  -  D  be  i  *  1, . . . ,  k  (thus  u\  =  z, 
and  follows  u,w,).  Since  P,  Q  Do  $  D, 

y(«.)  +  vK)  <  a  +  1, 

v(u>i)  +  y(m-H)  >  -l. 

Note  that  y(tii)  =  A;  y(u*)  <  6  (by  (4)  and  y{wk)  =  0);  and  y(u*+i)  >  y(uj  -  6  (also  by  (4)). 
These  three  inequalities  imply  that  for  any  j  €  Jl..fA/6  4 1]],  Px  has  some  w,  €  Uj.  For  a  given  j 
choose  the  last  such  t.  Then  u,+i  €  Uj- j.  Together  with  (4)  this  implies  w,  €  Wj,  since 

-bj  <  -yK+i)  -  b  <  y(w<)  -  o  -  1  <  -y(u,)  <  -b(j  -  1). 

We  have  shown  that  w,  €  W}  and  u,  €  Uj,  as  desired.  I 

Before  continuing  we  give  a  useful  refinement  of  Lemma  3.2.  Let  X  be  a  matching  such  that 
every  edge  not  in  X  has  multiplicity  at  most  Mx- 

Corollary  3.1.  For  some  constant  c,  at  any  point  in  match,  Av/$  <  eny/Mx- 

Proof.  The  proof  b  similar  to  the  lemma.  We  show  that  for  any  j  €  [l..f  A/6  4 1]],  each  P,  has 
an  edge  uw  not  in  X  with  u  €  Uj  U  Uj-% ,  w  6  Wj.  Thb  implies  Mx\Uj  U  Uj- 1  ||W/|  >  ♦,  which 
leads  to  the  desired  conclusion. 

To  find  the  desired  edge  uw  for  a  path  Pr,  proceed  exactly  as  in  the  lemma  to  find  an  index 
*  with  u,  6  Uj,  Ut+i  €  Uj- 1  and  w,  €  W}.  One  of  the  edges  u,w,,  w,u<+i  b  not  in  X  and  can  be 
taken  as  uw.  I 


Another  bound  on  A  is  useful  for  large  multiplicities.  It  is  similar  to  the  bound  used  in 
|GoT87a].  It  justifies  using  the  array  Q[l..dn]  to  compute  minima  in  the  Hungarian  search. 

Lemma  3.2.  For  some  constant  d,  at  any  point  in  matcht  A  <  in. 

Proof.  The  proof  of  Lemma  3.2  shows  that  for  any  €  [l..[A/6  +  1]],  P,  has  some  ti  €  Vj.  I 

Corollary  3.2.  The  number  of  iterations  of  match  is  0(min{\/F,n2^®A/1^8,n}). 

Proof.  Each  execution  of  Step  1  (except  possibly  the  first)  augments  along  at  least  one  path,  i.e., 
it  decreases  $  by  at  least  one.  The  definition  of  Step  1  implies  that  each  Hungarian  search  (except 
the  last)  increases  A  by  at  least  one.  Now  the  first  two  bounds  of  the  lemma  follow  because  at  any 

point  in  the  algorithm  A  or  ♦  is  at  most  B,  where  Lemma  3.1  gives  B  =  \/W  and  Lemma  3.2 

gives  B  -  (cn)2/3A/^a.  The  third  bound  follows  from  Lemma  3.3. 1 

The  corollary  implies  the  following  time  estimates.  The  estimates  are  good  for  graphs  or 
multigraphs  of  very  small  multiplicity. 

Theorem  3.1.  A  minimum  perfect  DCS  on  a  bipartite  multigraph  can  be  found  in 
0(min(v/^,«2^3M1^}7Jllog(nfV))  time.  The  space  is  O(m).  I 

For  example  in  a  bipartite  graph  a  minimum  perfect  DCS  can  be  found  in 
OfminfY/m.n^^TnlogfnlV))  time. 

The  bounds  of  the  theorem  also  apply  to  finding  a  minimum  cost  DCS.  To  see  this  let  G  be  the 
given  multigraph  or  graph.  Form  77  by  taking  two  copies  of  G  and  adding  a  set  of  edges  X,  where 
for  each  v  €  V'(G),  X  contains  an  edge  joining  the  two  copies  of  v,  with  multiplicity  u(v)  -  t(v)  and 
cost  tero.  It  is  easy  to  see  that  77  is  bipartite,  and  a  minimum  perfect  DCS  on  77  gives  a  minimum 
cost  DCS  on  G.  Furthermore  X  is  a  matching,  so  Corollary  3.1  applies  with  Mx  =  M.  This  implies 
the  time  bound  of  the  theorem  for  minimum  cost  DCS.  Minimum  cost  maximum  cardinality  DCS 
is  similar. 

Now  we  derive  bounds  that  are  good  for  multigraphs  with  moderately  sized  multiplicities. 
First  observe  Lemma  2.3  still  holds:  in  match  there  is  no  alternating  cycle  of  eligible  edges.  The 
proof  is  essentially  the  same:  There  is  no  such  cycle  initially,  since  the  edges  in  D  are  ineligible.  A 
Hungarian  search  does  not  create  such  a  cycle,  since  immediately  after  a  dual  adjustment  a  cycle 
leaving  7  on  a  new  eligible  edge  re-enters  7  on  an  ineligible  edge. 
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This  fact  ensures  that  the  time  for  a  depth-first  search  in  Step  1  is  0(m)  phis  the  totai 
augmenting  path  length.  Thus  the  total  time  for  match  is  0(mB  +  A ),  where  B  is  the  number  of 
iterations  and  A  is  the  total  augmenting  path  length.  Corollary  3.2  bounds  B\  now  we  estimate  A. 

Lemma  3.4.  A  =  0{mm{UloiU,nJSiU}). 

Proof.  As  in  Corollary  2.1,  A  <  U  +  A,.  For  the  first  bound  estimate  the  summation  as  in 

Corollary  2.1,  using  Lemma  3.1.  For  the  second  bound,  Lemma  3.2  shows  that  the  summation  is 
at  most  £f=1  cny/MJi  *  0(ny/MU).  I 

Theorem  3.2.  The  transportation  problem  (capacitated  or  not)  can  be  solved  in 
0((min{v^,n2/sJWJ/d,n}m  +  min{U\oiU,ny/Xflr))\o%(nN))  time.  The  space  is  0[m).  I 

To  understand  this  rather  involved  time  bound,  first  note  that  the  terms  containing  M  are 
relevant  only  in  the  capacitated  transportation  problem.  The  main  use  of  the  theorem  in  this 
paper  is  when  V  =  O(nm),  in  which  case  the  time  is  0(nmlogn)o%(nN));  this  bound  is  used 
in  Section  3.2  to  solve  transportation  problems  with  larger  U.  Fbr  further  applications  we  con¬ 
centrate  on  the  range  M  *  O(n).  In  this  case  the  above  bound  for  U  =  O(nm)  holds,  and  also 
the  bound  0(n2v/mlog(nAr));  hence  in  this  range  the  performance  is  competitive  with  [GoT87a]. 
In  most  of  the  range  M  =  O(n)  the  bounds  of  Theorem  3.2  are  those  of  Theorem  3.1  with  in 
replaced  by  m:  Using  U  log  U  as  the  second  term  of  the  time  bound  and  writing  Bm  as  the  first 
term,  the  first  term  dominates  if  V  «  0($m/ log  n).  Hence  the  bound  is  0(n2/8Ml/sm  log(n  A)) 
if  V  =  logn),  e.g.,  M  *  0(n/(logn)3/2);  the  bound  is  0(y/mMmlog(nN))  if 

V  =  0(iJmMm/\ogn),  e.g.,  M  ~  0(m/(logn)2);  the  bound  is  0( y/nMm log(nA^))  if  U  — 
0(\/nKimj logn),  e.g.,  all  degree  constraints  are  0{M)  and  M  =  0(m2/(n(logn)2). 

As  in  Theorem  3.1,  the  same  bounds  hold  for  networks  where  each  node  has  an  upper  and 
lower  bound  on  its  desired  degTee,  and  the  objective  is  minimum  coat  or  minimum  cost  maximum 
cardinality. 

3.2.  Scaling  edge  multiplicities. 

In  multigraphs  with  large  multiplicities  efficiency  is  gained  by  scaling  the  multiplicities.  Let 
D  be  a  DCS.  Recall  that  for  an  edge  t,  u(e)  and  d(e)  denote  the  multiplicities  of  e  in  G  and  D, 
respectively;  for  a  vertex  v,  u(v)  and  d(v)  denote  the  degree  constraint  of  t>  and  the  degree  of  v 
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in  D,  respectively.  The  term  u-value  refers  to  a  multiplicity  u(e)  or  a  degree  constraint  u(v).  The 
approach  is  to  scale  v- values.  The  ‘closeness  lemma*  needed  for  scaling  is  the  following. 

Let  G  be  a  multigraph  with  u-values  for  which  D  is  a  minimum  cost  maximum  cardinality 
DCS.  Form  u+  by  adding  one  to  the  u-values  of  an  arbitrary  subset  of  vertices  and  edges  (in 
particular  a  u-value  can  increase  from  zero  to  one).  Let  7  be  the  number  of  increased  u-values  (so 
/  <  m  +  n).  Let  D+  be  a  minimum  cost  maximum  cardinality  DCS  for  u+.  Let  77+  ©  D  denote 
the  subgraph  that  is  the  direct  sum  of  subgraphs  D+  and  D  (i.e.,  for  any  edge  e,  D+  ©  D  has 
|d+(e)  -  d(e)|  copies  of  e).  Choose  D*  so  that  \D+  ©  D\  is  as  small  as  possible.  ^+(e,  D)  denotes 
the  deficiency  of  D  at  v  for  u-values  u+. 

Lemma  3.6.  Z?+  ©  D  can  be  partitioned  into  at  most  I  simple  alternating  paths  and  cycles 
(where  ‘alternating*  means  edges  are  alternately  in  D  and  D+). 

Proof.  Since  both  D+  and  D  are  DCS’s  for  u+,  D+&D  can  be  partitioned  into  simple  alternating 
paths  and  cycles;  for  each  vertex  v,  at  most  4>+(v,  Z?+)  paths  end  at  v  on  a  D-edge,  and  similarly 
for  a  TH-edge.  Call  an  edge  vu>  with  d+(vw)  >  d(vw)  new  if  either 
(i)  d(vto)  ss  v (vw),  or 

(it)  vw  is  an  end  of  a  path  of  D+  ©  D  and  d(v)  =  u(v). 

There  are  at  most  I  new  edges.  (A  type  (t)  new  edge  clearly  has  an  increased  u-value.  For  a  type 
(»i)  new  edge  vw,  v  has  an  increased  u-value  and  ^+(v,D)  «  1,  so  vw  is  the  only  type  (it)  edge 
associated  with  v).  Thus  it  suffices  to  show  that  any  alternating  path  or  cycle  P  of  D+  ©D  contains 
a  new  edge. 

P  does  not  begin  and  end  with  a  D-edge,  since  D+  has  maximum  cardinality.  Suppose  P  does 
not  contain  a  new  edge.  Then  D  ©  P  is  a  feasible  DCS  for  u.  (This  follows  since  a  D+-edge  vw  of 
P  has  d(vu/)  <  u(vw)\  further  if  this  edge  vu>  is  an  end  of  P  then  d(v)  <  u(v)).  P  does  not  begin 
and  end  with  a  D+-edge,  since  D  has  maximum  cardinality.  Thus  P  is  an  even  length  alternating 
path  or  cycle.  This  implies  P  has  tero  net  cost  (with  respect  to  D  or  D+).  But  this  contradicts 
the  fact  that  \D+  ©  D\  is  as  small  as  possible.  I 

The  lemma  indicates  that  D+  ©  D  can  be  found  in  a  'small*  multigraph  G',  defined  as  follows. 
A  vertex  v  €  V(C)  corresponds  to  vj,  V2  €  V(Gf);  G'  has  an  edge  v\vz  of  cost  0  and  multiplicity 
I.  An  edge  vw  6  E(G)  corresponds  to  edges  vjwj,  v^v>2  £  E(G'),  with  multiplicities  and  costs 

u'(viu>i)  *  u+(tw)  -  d(vw),  e'(viu>i)  =  c[vw)  -  nN ; 
u'(v2U>2)  =  d(vw),  c'(v2W2)  =  -e[vw)  +  nN. 


Call  these  edges  G-edges\  edges  v?w?  are  D-edges  and  edges  vjtt>i  are  non- D- edges.  Finally  each  v  6 
V(G)  has  upper  and  lower  degree  constraints  u'(vi)  =  u'(v2)  =  /,  =  0,  t (V2)  =  7-^+(t>,  D). 

Consider  D a  minimum  cost  DCS  of  G '.  It  is  easy  to  see  that  the  G -edges  of  D'  can  be 
partitioned  into  at  most  7  paths  and  cycles  that  are  alternating  for  D,  and  that  D©  D'  is  a  feasible 
DCS  (note  that  the  lower  bounds  in  G'  allow  a  vertex  t>  to  be  on  at  most  4>+(t\  D)  non- D-edges  at 
the  end  of  an  alternating  path).  Furthermore,  the  costs  of  G -edges  guarantee  that  D’  has  as  many 
augmenting  paths  as  possible  and  no  *de- augmenting  paths*  (i.e.,  paths  that  begin  and  end  with 
D-edges).  Thus  D'  is  the  desired  set  of  paths  (D+  ©  D  or  an  equivalent  set). 

Now  we  state  the  capacity  scaling  algorithm  for  finding  a  minimum  perfect  DCS.  Given  a  DCS 
problem  on  a  multigraph  G,  let  u  denote  the  given  v-values,  with  M  the  largest  «-value.  (Without 
loss  of  generality  M  is  the  V-value  of  a  vertex).  Consider  each  V- value  to  be  a  binary  number 
61 . . .  bk  of  it  =  [log  MJ  +  1  bits.  The  routine  maintains  u  as  the  ti-values  in  the  current  scale,  and 
D  (and  d)  as  the  DCS  in  the  current  scale.  It  initialises  each  u(e),  d(e)  and  u(v)  to  zero.  Then  it 
executes  the  following  loop  for  scale  index  s  going  from  1  to  k: 

Step  1.  For  each  edge  e,  d(c)  «-  2d(e)  and  u{c)  —  2u(e)-f  (bit  b ,  of  0(e)).  For  each  vertex 
v,  u(v) «-  2u(v)+  (bit  b,  of  ff(t/)). 

Step  2.  Form  the  multigraph  G'  described  above.  (Note  1  <  m  4-  n;  u+  is  given  by  u  in  the 
algorithm.) 

Step  3.  Let  D'  be  a  minimum  cost  DCS  on  G'.  Set  D  —  D  ©  D'  and  let  d  be  the  function 
corresponding  to  D.  I 

The  correctness  of  this  algorithm  follows  from  the  above  discussion.  (Note  that  this  algorithm 
works  on  both  bipartite  and  general  graphs).  To  analyze  the  running  time,  assume  that  the  DCS  in 
Step  3  is  found  using  the  algorithm  of  Theorem  3.2  (note  that  V  —  O(nm)).  The  following  bound 
is  immediate. 

Theorem  3.3.  The  transportation  problem  (capacitated  or  not)  can  be  solved  in 
0(nm  log  n  log(n^)  log  M)  time.  The  space  is  0(m).  I 

This  result  extends  to  the  variants  of  the  perfect  DCS  problem  mentioned  above. 

The  capacity  scaling  algorithm  has  a  variant,  called  EK  capacity  scaling,  since  it  is  similar 
to  the  mincost  flow  algorithm  of  Edmonds  and  Karp  (EK).  EK  capacity  scaling  is  used  below  to 
solve  flow  problems  with  lower  bounds.  It  maintains  a  set  of  dominated  tight  duals  on  G  for  v  and 


D.  (“Dominated  tight*  duals  satisfy  the  inequalities  for  l-feasiblity  with  cost-length  cl  replaced 
by  cost  c).  The  algorithm  scales  up  costs  and  duals  exactly  as  in  Step  1  of  the  main  routine  and 
scalejmatch.  This  ensures  that  the  costs  input  to  match  are  positive  integers  except  for  edges  that 
were  in  the  DCS  of  the  previous  scale. 

The  match  routine  initializes  all  duals  y{v)  to  0  and  the  DCS  to  (e|c(e)  <  0}.  This  DCS  does 
not  violate  any  degree  constraint  and  the  duals  are  dominated  tight.  Then  match  does  as  many 
minimum  cost  augmentations  as  possible.  The  augmenting  paths  are  found  using  the  Hungarian 
algorithm.  Then  match  searches  for  alternating  even  paths  and  cycles  with  minimum  negative 
cost,  and  augments  along  them  (“augmenting  along  a  path  or  cycle”  is  defined  by  analogy  with 
augmenting  along  an  augmenting  path;  such  an  augmentation,  in  the  current  context,  gives  a  new 
matching  of  equal  cardinality).  These  searches  are  also  implemented  with  the  Hungarian  algorithm. 
(A  minor  point  is  that  match  uses  the  standard  Hungarian  algorithm,  based  on  dominated  tight 
duals  as  in  [L].  This  differs  from  the  Hungarian  search  in  Section  3.1  which  uses  1-feasible  duals.  The 
main  difference  is  the  definition  of  “eligible”,  which  uses  either  cost  or  cost-length  as  appropriate). 

Each  Hungarian  search  takes  time  0(m  +  n  log  n)  using  Fibonacci  heaps  [FT].  In  problems 
where  each  scale  has  I  -  O(n),  the  total  time  for  EK  capacity  scaling  is  0(n(m  +  n logo) log  Af), 
slightly  improving  Theorem  3.3.  (The  correctness  of  EK  capacity  scaling  follows  from  Lemma  3.5, 
applied  one  edge  or  vertex  at  a  time,  i.e.,  I  =  1). 

Next  consider  the  transportation  problem  urith  cost  functions.  This  problem  allows  parallel 
copies  of  an  edge  to  have  different  costs.  Specifically  the  cost  of  the  p‘*  copy  of  an  edge  e,  1  < 
V  -  v(e)i  “  t>iveD  ty  c(e,p),  a  nondecreasing  function  of  p  that  can  be  evaluated  in  0(1)  time.  Ab 
usual  these  costs  are  in  [-.N.../V],  and  each  vertex  v  has  a  desired  degree  u(i>).  The  problem  is  to 
find  a  minimum  cost  perfect  DCS  for  these  degree  constraints.  Note  that  the  desired  DCS  can  still 
be  represented  by  an  integral  function  on  the  edges  d{c),  where  0  <  d(c)  <  u(e),  since  without  loss 
of  generality  the  DCS  contains  the  d(c)  copies  of  c  with  smallest  cost. 

As  examples  of  this  problem,  c(e,p)  =  [a«pj  +  6«  is  the  original  DCS  problem  for  a,  =  0  and 
a  simple  example  of  diminishing  returns  to  scale  for  ae  >  0.  Alternatively  r(e,p)  could  be,  say,  a 
piecewise  quadratic  function;  in  this  case  evaluating  c(e,p)  for  arbitrary  p  would  probably  involve 
a  binary  search  on  the  breakpoints.  (Note  that  in  the  definition  of  the  transportation  problem 
with  cost  functions,  the  restriction  to  nondecreasing  cost  functions  e(e,p)  is  crucial:  without  it  the 
problem  is  NP-hard  [GJ,  p.214].  Also  note  that  the  solution  to  the  problem  is  a  multigraph  with 
integral  multiplicities,  by  definition.  This  assumption  of  integrality  is  also  crucial.  This  issue  is 
discussed  further  after  Theorem  3.5  below  for  network  flows,  where  real-valued  flows  make  sense). 
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The  algorithm  of  Theorem  3.3  works  correctly  for  the  transportation  problem  with  cost  func¬ 
tions,  if  we  treat  parallel  copies  of  an  edge  with  different  costs  as  different  edges.  In  the  time  bound, 
however,  the  term  m  now  counts  each  edge  e  according  to  the  number  of  distinct  costs  c(e,p).  We 
will  show  how  to  extend  the  capacity  scaling  algorithm  to  the  transportation  problem  with  cost 
functions,  preserving  the  time  bound  of  Theorem  3.3. 

First  we  modify  the  cost  scaling  algorithm  to  preserve  the  time  bounds  of  Theorems  3.1-2.  The 
derivation  of  those  theorems  remains  valid  for  cost  functions  and  gives  the  desired  time  bounds, 
provided  all  individual  steps  are  implemented  to  run  in  essentially  the  same  time  as  before.  This 
means  implementing  Step  1  of  the  main  routine  and  scale-match  in  time  O(m)  (even  though  they 
modify  every  cost)  and  similarly  for  match.  This  can  be  done  because  of  the  following  observation. 
The  conditions  for  a  DCS  D  to  be  1-feasibile  are  equivalent  to  a  system  with  only  two  inequalities 
per  edge  e  =  vw, 

c(e,d(e))  <  y[v)  +  y(w)  <  c(e,d(e)  +  1)  +  1. 

Further  the  only  copies  of  e  that  can  be  eligible  are  D -edges  costing  c(c,d(e))  and  non- Z) -edges 
costing  c{e,d(e)  +  1)  +  1. 

Step  1  of  the  main  routine  and  scale-match  do  not  explicitly  modify  edge  costs.  Instead  match 
computes  the  cost  of  an  edge  when  it  is  needed  in  0(1)  time  using  arithmetic.  Specifically,  the  ptk 
cost  for  vw  is 

tr*nc(r(vu,  p)/2*-')  -  yo(v)  -  yoM,  (5) 

where  tnme  denotes  integer  truncation,  k  =  [log(n  +  2)tfJ  is  the  number  of  cost  scales,  s  is  the 
index  of  the  current  cost  scale,  and  yo  denotes  duals  at  the  start  of  scale  s. 

The  match  routine  starts  by  initializing  D  to  contain  all  edges  costing  less  than  -1.  This  is 
done  by  examining  each  edge  and  adding  smallest  cost  copies  to  D  until  the  cost  reaches  -1.  The 
time  is  0(m  +  U),  which  suffices  for  the  bounds  of  Theorems  3.1-2. 

In  the  depth-first  search  of  Step  1,  it  is  unnecessary  to  know  the  multiplicity  of  each  eligible 
edge  when  the  search  begins.  Rather,  costs  e(e,  d(e))  and  c(e,  d(e)  + 1)  are  used  to  determine  which 
edges  have  at  least  one  eligible  copy.  When  the  depth-first  search  finds  an  augmenting  path  P,  the 
next  cost  for  each  edge  e  6  P  is  used  to  see  if  there  is  another  eligible  copy  of  e  (i.e.,  for  e  £  PHD, 
another  copy  of  e  is  eligible  if  e(e,d(e)  - 1)  =  c(e,d(e)),  and  similarly  for  e  €  P-  D).  Thus  the  time 
for  the  depth-first  search  is  still  0(m)  plus  the  total  augmenting  path  length.  It  is  obvious  that 
the  Hungarian  search,  given  costs  e(e,d(e))  and  e(e,d(e)  +  1),  uses  time  0(m).  Thus  the  bounds 
of  Theorems  3.1-2  apply. 


Now  we  modify  the  capacity  scaling  algorithm  of  Theorem  3.3.  The  new  version  works  by 
scaling  the  domain  of  the  cost  functions.  The  closeness  lemma  (Lemma  3.5)  generalises  as  follows. 
Let  G  be  a  multigraph  with  c  functions  and  u  values  for  which  D  is  a  minimum  cost  maximum 
cardinality  DCS.  Form  u+  by  adding  one  to  the  u-values  of  an  arbitrary  set  of  vertices  and  edges. 
Form  c+  so  that  for  each  edge  e  and  p  £  (0..u(e)), 

e(e,p+l)>e+(c,p  +  l)>e(c,p)-  (6) 

(Here  c(e,0)  =  -oo).  Let  /  be  the  number  of  vertices  with  an  increased  u- value  plus  the  number 
of  edges  with  an  increased  u- value  or  some  decreased  e- value  (so  that  1  <  m  +  n).  Let  £>+  be  a 
minimum  cost  maximum  cardinality  DCS  for  c+  and  ti+,  chosen  so  that  \D+  ©  D\  is  minimum 
(D+  ®  D  has  the  obvious  interpretation). 

Lemma  3.6.  D+  ©  D  can  be  partitioned  into  at  most  I  simple  alternating  paths  and  cycles. 

Proof.  The  argument  is  an  expanded  version  of  Lemma  3.5;  we  will  give  only  the  new  material. 
The  definition  of  new  edge  is  expanded  to  include  a  type  (m)  new  edge  e,  defined  to  have  rf+fe)  > 
d(e)  and  c(e,d(e)  +  1)  >  c+(e,d(c)  +  1),  where  only  the  d(e)  +  1“  copy  of  e  is  new.  (Note  that 
cf+(e)  -  d(e)  may  be  larger  than  one). 

The  argument  expands  in  the  case  that  P  is  an  even  length  alternating  path  or  cycle  not 
containing  a  new  edge,  and  we  must  show  that  it  has  zero  net  cost  (with  respect  to  c+  and  D+). 
The  net  cost  of  P  with  respect  to  c+  and  D+  is  nonnegative,  by  the  choice  of  D+.  Hence  it  suffices 
to  show  that  the  net  cost  of  P  with  respect  to  e+  and  D  is  nonnegative. 

This  follows  from  the  choice  of  D  if  for  every  edge  e  whose  plh  copy  is  in  P  D  D+,  c+(e,p)  > 
c(e,  <f(e)  +  1).  We  prove  this  inequality  as  follows.  The  copy  of  e  is  not  new  and  p  >  d(e)  +  1.  Now 
consider  two  cases:  If  p  =  d{e)  +  1  then  c+(e,p)  =  c(e,d(e)  +  1),  as  desired.  If  p  >  d[e)  +  1  then 
c+(«»p)  >  c(e,d(e)  +  1)  by  (6),  as  desired.  I 

This  lemma  justifies  an  algorithm  similar  to  the  capacity  scaling  algorithm.  The  main  dif¬ 
ferences  are  as  follows.  Step  1,  in  addition  to  scaling  d  and  u,  scales  the  cost  function  domain. 
Specifically  let  co  denote  the  given  cost  function.  Then  for  each  t  and  p  £  [l..u(e)]  (where  u(e) 
is  the  new  u-vahie)  Step  1  sets  e(e,p)  •-  co(e,  (2fc-*pJ),  where  k  =  [log Afj  +  1  is  the  number  of 
capacity  scales  and  »  is  the  index  of  the  current  capacity  scale.  Observe  that  the  DCS  for  the 
new  d  is  a  minimum  cost  maximum  cardinality  DCS  for  (the  new)  u-values  rounded  down  to  even 
numbers  and  (the  new)  c(e,2p  -  1)  increased  to  c(e,2p).  So  Lemma  3.6  applies  and  justifies  the 
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remaining  steps:  Step  2  defines  G'  as  before  but  with  costs  changed  in  the  obvious  way  to  take  cost 
functions  into  account.  Step  3  computes  D'  using  the  cost  scaling  algorithm  described  above. 

For  efficiency,  these  three  steps  are  not  done  explicitly.  (For  instance,  doing  Step  2  explicitly 
would  use  9(m2)  time,  since  an  edge  can  be  in  G'  with  multiplicity  m  +  n).  Step  1  computes  only 
two  new  costs  for  each  edge  e,  c(c,d{c))  (already  known)  and  e(e,d(e)  +  1).  To  do  Step  2,  G'  is 
initialized  to  contain  only  the  cheapest  copy  of  each  edge  of  type  viwi,  V2tt^.  This  is  the  copy  that 
will  be  added  to  the  DCS  D’  first.  Each  copy  comes  from  a  cost  computed  in  Step  1.  When  the 
cost  scaling  algorithm  checks  to  see  if  there  is  another  eligible  copy  of  an  edge  e  (in  the  depth-first 
search),  the  next  higher  (or  lower)  cost  copy  of  t  is  computed,  its  cost  is  scaled  down  using  (5)  and 
it  is  used  in  the  cost  scaling  algorithm. 

Theorem  3.4.  The  transportation  problem  (capacitated  or  not)  with  cost  functions  can  be  solved 
in  0(nm  log  nlog(nlV)  log  M)  time.  The  space  is  0(m).  I 

3.3.  Network  flow. 

Our  results  extend  to  integral  network  flows.  It  is  convenient  to  work  with  the  problem 
of  finding  a  minimum  cost  circulation,  defined  as  follows  |L].  Let  G  be  a  directed  graph  with 
nonnegative  integral  capacities  u(v)  for  each  vertex  v,  and  for  each  edge  e,  nonnegative  integral 
capacities  u(e),  lower  bounds  i(t)  and  costs  e(e).  The  minimum  circulation  problem  is  to  find  a 
feasible  circulation  with  smallest  possible  cost.  (If  vertex  capacities  are  not  given,  setting  u(v)  = 
«(vty)  does  not  change  the  problem.  The  circulation  problem  includes  the  minimum  cost  flow 
problem  as  a  special  case.  As  already  mentioned,  the  usual  definition  of  the  circulation  (network 
flow)  problem  allows  real- valued  parameters.  However  note  that  if  all  capacities  and  lower  bounds 
are  integral,  an  optimum  circulation  (flow)  that  is  integer-valued  always  exists  [L]). 

A  minimum  circulation  problem  on  a  network  G  can  be  transformed  to  a  minimum  perfect  DCS 
problem  on  a  bipartite  multigraph  B,  as  follows.  A  vertex  v  €  V(G)  corresponds  to  vj,  V2  6  V(£); 
B  has  an  edge  V]«2  of  cost  0  and  multiplicity  u(v).  An  edge  w  €  E(G)  corresponds  to  vi  u>2  €  E{B) 
with  cost  c(w)  and  multiplicity  u(vtc)  -  l(vw).  The  degree  constraints  on  B  are 

u(vi)  =  «(t/)  -  ^{^vw)|tw  €  E(G)}, 
v(vz)  =  u(v)  -  J"^{l(igi>)[tPi?  €  E(G)}. 

A  circulation  on  G  corresponds  to  a  perfect  DCS  on  B  costing  less  by  exactly 
£{£(e)e(e)|e  €  E[G)).  Thus  the  flow  problem  can  be  solved  using  the  DCS  algorithms  given 


above.  Note  that  B  has  n  vertices  in  each  vertex  set,  O(m)  edges,  U  =  0(£{«(v)lv  €  V'(G)}) 
and  frt  »  0(U  +  £{«(e)|e  €  £(£)}.  In  part  (a)  below,  m  is  the  number  of  edges,  with  each  edge 
coanted  according  to  its  capacity. 

Theorem  3.5.  A  minimum  cost  circulation  on  a  network  with  all  edge  capacities  and  lower 

bounds  in  [0 _ A/]  can  be  found  in  the  following  time  bounds  (and  space  O(m)): 

(o)  0(min{v/W,n2^8A/1^8}ffflog(nlV)). 

(6)  ,n)m  +  min{mA/log(mA/),n2 y/rn})\oz{nX)). 

(c)  0(nm  log  n)og(nAf)  log  M). 

These  bounds  also  hold  when  each  edge  cost  is  a  convex  function  of  its  flow. 

Proof.  These  bounds  follow  essentially  from  Theorems  3.1*4.  Note  that  M  does  not  necessarily 
bound  the  multiplicities  in  B,  since  we  assume  no  bound  on  vertex  capacities  in  G.  Nonetheless  the 
bound  for  part  (b)  holds:  To  show  this  use  Corolla^  3.1,  with  matching  X  containing  all  edges  of 
the  form  note  that  Mx  -  M.  Also  the  bound  for  part  ( e)  holds:  There  are  log(mM)  capacity 
scales,  but  the  time  bound  involves  the  factor  log  M  because  the  first  log  m  scales  are  trivial.  I 

Note  that  in  Theorem  3.5  (c),  the  algorithm  for  convex  cost  functions  finds  an  optimal  integral* 
valued  flow.  However  this  flow  need  not  be  the  global  optimum,  which  may  involve  real* valued  flow 
values.  Finding  this  solution  appears  to  be  much  harder.  For  instance  if  the  cost  of  an  edge  is  a 
quadratic  function  of  its  flow,  finding  a  minimum  cost  flow  is  NP*hard  [GJ,  H]. 

Next  consider  a  minimum  circulation  problem  in  which  0(n)  edges  have  finite  capacity  (every 
edge  still  has  a  lower  bound,  perhaps  zero).  Such  problems  arise  as  covering  problems;  a  common 
special  case  is  circulations  with  lower  bounds  but  no  upper  bounds  (e.g.,  the  aircraft  scheduling 
problem  of  (L,  p-139]). 

Theorem  3.6.  A  minimum  circulation  on  a  network  with  lower  bounds  but  only  O(n)  finite 
capacities,  all  in  (O..A/j,  can  be  found  in  0(n(m  +  n  log  n)  log(nAf))  time  and  0{m)  space. 

Proof.  Without  loss  of  generality  assume  that  no  cycle  has  negative  cost  and  infinite  capacity. 
Then  it  is  easy  to  see  that  all  infinite  capacities  (on  edges  or  vertices)  can  be  replaced  by  any 
number  that  is  at  least  5  *  £{tf  e(e)  is  finite  then  c(e)  else  t(e)\e  €  £((?)}  +  £{c(v)|t>  € 
V(G),  e(v)  is  finite}. 

The  algorithm  is  as  follows:  Find  5  and  k  =  flog  5].  For  each  infinite  capacity  vertex  v,  set  its 
capacity  to  5;  for  each  infinite  capacity  edge  e,  set  its  capacity  to  £(e)  +  2*+1 .  Transform  the  new 


circulation  problem  to  a  DCS  problem,  as  above,  and  solve  the  DCS  problem  using  EK  capacity 
scaling. 

The  correctness  of  this  algorithm  follows  from  the  opening  remark.  To  estimate  the  efficiency, 
note  that  in  the  DCS  problem,  every  infinite  capacity  edge  of  G  has  multiplicity  2*+’  and  every 
vertex  vi,  has  degree  constraint  at  most  S  <  2*.  Hence  the  first  scale  is  trivial—  no  edges  are 
in  the  DCS,  and  the  duals  can  be  set  to  any  values  small  enough  so  that  they  are  dominated  on  all 
edges.  Every  scale  after  the  first  has  I  =  0(n)  (recall  that  I  is  the  number  of  increased  upvalues), 
since  the  upvalues  of  infinite  capacity  edges  double.  Hence  the  total  time  is  0{n(m  +  n  log  n)  log  S), 
implying  the  desired  bound.  I 

Corollary  S.3.  The  directed  Chinese  postman  problem  can  be  solved  in  0(nm  log  n)  time  and 
O(m)  space. 

Proof.  We  use  the  notation  for  the  postman  problem  given  in  the  Introduction.  Transform  the 
postman  problem  to  a  circulation  problem  as  follows:  Add  new  vertices  s  and  t.  For  each  start 
vertex  v  add  edge  av  with  capacity  indegree(v)  -  outdegree(v);  for  each  end  v  add  edge  vt  with 
capacity  outdegree(v)  -  indegree[v).  If  q  is  the  value  of  a  maximum  flow  in  this  graph,  add  edge  ta 
with  lower  bound  q.  (9  is  the  cardinality  of  the  desired  set  of  paths  P).  Define  the  capacity  of  each 
vertex  and  each  original  edge  of  G  as  infinite.  The  problem  is  to  find  a  minimum  cost  circulation 
from  s  to  t. 

The  correctness  of  this  transformation  is  obvious.  To  estimate  the  efficiency,  first  observe  that 
q  can  be  found  in  O(nm)  time  (in  any  network  with  integer  capacities  and  maximum  flow  value 
O(m),  a  maximum  flow  can  be  found  in  O(nm)  time  [G85a]). 

The  circulation  problem  can  be  solved  in  time  0(n(m+n  log  n)  log  n),  by  Theorem  3.6.  For  m  > 
n  logn  this  is  the  desired  bound.  If  m  <  n  logn  the  desired  time  is  achieved  by  any  algorithm  that 
does  O(m)  Hungarian  searches.  For  instance  the  desired  circulation  can  be  found  without  scaling, 
by  repeatedly  finding  a  minimum  cost  augmenting  path  from  a  to  t  (the  number  of  augmentations 
is  at  most  m).  Alternatively  the  DCS  graph  for  Theorem  3.6  can  be  modified  slightly  so  that 
the  scaling  algorithm  does  at  most  m  augmentations,  thereby  achieving  the  desired  bound  in 
this  case  too.  (The  idea  is  to  define  capacities  so  that  the  scaling  algorithm  does  at  most  one 
search  per  new  unit  of  vertex  capacity.  Specifically,  for  a  start  vertex  v,  in  the  DCS  graph  define 
u(vi)  =  2*  +  tndegree(v)  -  o\itdtgret{v),  u(vj)  =  2*;  similarly  for  an  end  vertex  or  other  vertex). 
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4.  Conclusions. 

Table  1  shows  that  in  terms  of  asymptotic  estimates,  many  network  problems  can  be  solved 
efficiently  by  scaling.  Scaling  algorithms  also  tend  to  be  simple  to  program.  For  instance  the 
assignment  algorithm  consists  of  an  outer  scaling  loop  plus  an  inner  loop  that  does  a  depth-first 
search  followed  by  a  Dijkstra  calculation.  We  believe  that  such  algorithms  will  run  efficiently  in 
practice.  Note  that  in  the  experiments  done  by  Bateson  [Ba]  the  scaling  algorithm  of  [G85a]  beat 
the  Hungarian  algorithm  as  long  as  the  cost  of  the  matching  could  be  stored  in  a  machine  integer. 
Our  assignment  algorithm  has  even  simpler  code  than  [G8Saj  and  so  should  do  even  better. 

The  assignment  algorithm  has  a  processor-efficient  parallel  implementation.  Details  are  given 
in  (GabT87b].  Also,  the  assignment  algorithm  extends  to  matching  on  general  graphs.  The  time 
bound  for  finding  a  minimum  perfect  matching  on  a  general  graph  is  0(yno(m,n)lognmlog(nfV)). 
The  algorithm  is  more  complicated  because  of  ‘blossoms*  that  occur  in  general  matching,  which 
compound  the  error  due  to  scaling.  Details  are  in  (GabT87aj. 
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i: 


Strong  Polynomial  Bound 

Assignment  problem 
0(n(m  +  n  log  «))  [FT] 


New  Scaling  Bound 


Of^mloglnA/)) 


Shortest  paths  (single-source,  directed  graph,  possibly  negative  lengths) 

O(nm)  (Bel)  0(y/nm  log(titf)) 

Minimum  cost  degree-constrained  subgraph  of  a  bipartite  multigraph 


0(1/ (tn  +  n  log  n))  [FT,G83] 

Transportation  problem  (uncapacitated  or  capacitated) 
0(min{t/,nlogf/}(fn -fnlogn))  (FT,  EK,  L] 

Minimum  cost  flow 

0(n2(m  +  b  log  n)  log  n)  jGalT} 


Directed  Chinese  postman  problem  (with  edge  lengths) 
0(m(m  +  nlogn))  [FT,  EJ,  Le] 


0(min{y/U,  rt2/3AP/3}mlog(nN)) 


0((mm{y/P,n}m  +  l/loglf)log(nN)) 


0(nm  log  n  log(n  JV )  log  M) 
convex  cost  functions  allowed 
0(ft(m  +  n  log  n)  log  (nA/)) 
lower  bounds  only 


O(nmlogft) 


Table  I.  Bounds  for  network  problems. 

Parameters: 

n  «*  number  of  vertices  V  =  total  degree  constraints 

m  =  number  of  edges  N  =  maximum  cost  magnitude 

m  -  number  of  edges  counting  multiplicities  A/  =  maximum  flow  capacity  or  lower  bound, 

or  edge  multiplicity 
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